#pragma once
extern "C" {
#include <f2c.h>
#include <clapack.h>
}

namespace ClapackNet {
	using namespace System;

	public ref class DoubleGeneralMatrix {
	private:
		double * a;
		long * ipiv;
		int m_rownum;
		int m_colnum;

	protected:
		!DoubleGeneralMatrix() {
			delete[] a;
			delete[] ipiv;
		};
	public:
		DoubleGeneralMatrix(int rownum, int colnum) {
			m_rownum = rownum;
			m_colnum = colnum;
			a = new double[rownum * colnum];

			for (int i = 0; i < rownum * colnum; i++) {
				a[i] = 0.0;
			}
			ipiv = new long[rownum];
		};

		~DoubleGeneralMatrix(){
			this->!DoubleGeneralMatrix();
		};

		property int RowNum {
			int get() {
				return m_rownum;
			}
		};

		property int ColNum {
			int get() {
				return m_colnum;
			}
		};

		property double default[int, int] {
			double get(int i, int j) {
				return a[j * m_rownum + i];
			}
			void set(int i, int j, double value) {
				a[j * m_rownum + i] = value;
			}
		};

		void Inverse() {

			long n = m_rownum;
			long lda = n;
			long lwork = n;

			double * work = new double[lwork];
			long info = 0;

			dgetri_(&n, a, &lda, ipiv, work, &lwork, &info);

			delete[] work;

		};

		void Factorize() {

			long m = m_rownum;
			long n = m_colnum;
			long lda = m;
			long info;

			dgetrf_(&m, &n, a, &lda, ipiv, &info);
		};

		void Solve(DoubleGeneralMatrix ^ b) {

			char trans = 'N';
			long n = m_rownum;
			long nrhs = b->m_colnum;
			long lda = n;
			long ldb = n;
			long info = 0;

			dgetrs_(&trans, &n, &nrhs, a, &lda, ipiv, b->a, &lda, &info);
		};

		void FactorizeAndSolve(DoubleGeneralMatrix ^ b) {

			long n = m_rownum;
			long nrhs = b->m_colnum;
			long lda = n;
			long ldb = n;
			long info = 0;

			dgesv_(&n, &nrhs, a, &lda, ipiv, b->a, &ldb, &info);

		};
	};
}